Approximate, analytic solutions of the Bethe equation for charged particle range 
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By either performing a Taylor expansion or making a polynomial approximation, the Bethe equa- 
tion for charged particle stopping power in matter can be integrated analytically to obtain the range 
of charged particles in the continuous deceleration approximation. Ranges match reference data to 
the expected accuracy of the Bethe model. In the non-relativistic limit, the energy deposition 
rate was also found analytically. The analytic relations can be used to complement and validate 
numerical solutions including more detailed physics. 
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Keywords: charged particle, energy loss, stopping power, ion implantation 



INTRODUCTION 

The deceleration of ions as they pass through matter 
is important in a wide range of fields: medical ion radia- 
tion therapy, such as the treatment of tumors [![; radio- 
raphy with ions Q; radiolysis of chemical compounds 
; ion implantation in material processing and semi- 
conductor doping gas discharge plasmas locality 
of energy deposition in nuclear fusion plasmas, includ- 
ing ion beam heating for controlled thermonuclear fusion 
6 1 ; the design of radiation shielding for nuclear reactors 
7 j and spacecraft and particle physics experiments 
9|. Energy loss for different combinations of ion specie, 
ion energy, and decelerating material has been measured 
since the early 20th century [l(J, with a corresponding 
development in theoretical work. Calculations of ion en- 
ergy loss, and hence range, are very frequently performed 
using variants of the Bethe 11 1 and Bethe-Bloch [2] re- 
lations. Experimental and theoretical developments have 
focused on making corrections to the original Bethe re- 
lation to account for details of interactions with bound 



electrons and crystal structures at low energies [13j, [lj] , 
and quantum-mechanical limits to the transfer of energy 
under extreme relativistic conditions [HI]. The original 
Bethe relation is still used widely, particularly when a 
reliable, approximate result is needed rapidly, or in the 
fairly wide range of energies where the Bethe relation is 
adequately accurate [l^ ]. 

The Bethe relation describes the stopping power: the 
rate at which a moving ion loses energy to the surround- 
ing material. It is not trivial to use this relation to ob- 
tain the range of an ion, i.e. the distance for it to lose 
all of its kinetic energy. In practice, ion ranges are calcu- 
lated using numerical integration in multi-physics com- 
puter programs, or from scaling laws normalized to the 
range for other energies or masses. However, reliance on 
sophisticated computer programs for infrequent calcula- 
tions, without the ability to make a compact analytic 
estimate, can lead to errors. Analytic solutions are also 
valuable for validating computer programs reproducing 



the same physics. Here we point out analytic solutions 
to accurate approximations of the Bethe relation, which 
can give good estimates of ion ranges in matter. 



RANGE FROM STOPPING POWER 

In the continuous deceleration approximation, charged 
particles traversing matter lose kinetic energy E at a rate 
depending on their instantaneous energy and the local 
material. Expressed as the energy loss rate per distance 
traveled, the stopping power dE/dx can be used to de- 
termine the range I of the particle, by integrating the 
deposition until the particle is stationary: 
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However, dE/dx is expressed naturally in terms of E 
rather than x. Rearranging, 
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The integral can be found numerically for arbitrary stop- 
ping powers, or analytically for stopping powers of suffi- 
ciently simple form. 

The Bethe equation [ll[ describes the deceleration of 
charged particles by interaction with the electrons in 
matter: 



dE 
dx 



4tt NZz 2 ( q 2 



s c 2 (3 2 V47re 



111 



2m e c 2 p 2 



P 2 



(3) 

where (3 = v/c, v is the ion speed, I is the effective 
ionization of the target material, Z and z are the atomic 
numbers of the target and ion species respectively, N is 
the number density of target nuclei, m e is the mass of an 
electron, eo is the permittivity of free space, and c is the 
speed of light. 

To find the range of charged particles from Eq. [3l Eq. [2] 
can be integrated numerically, though this is not straight- 
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forward because of a singularity at low energies. We have 
not found an analytic solution for the integral. 



TAYLOR EXPANSION 

—dx/dE can be expanded as a Taylor series to make 
it more tractable for integration. This can be done for 
Eq.[3]with a relativistic expression for (3(E); we do it also 
for a non-relativistic (3(E) because the resulting integral 
is more amenable to subsequent manipulation. 



Relativistic 

The relativistic relation between (3 and kinetic energy 
E is 



(3(E) = 



y E(E + 2m,c 2 ) 
E + mid 2 



(4) 



where rrii is the rest-mass of the moving ion. For later 
convenience, we scale key quantities to be dimensionless. 
Substituting into Eq. [3] and expanding about zero, 
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Integrating Eq. the range is 
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where Ei(z) is the exponential integral function, 



dt. 
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Non-relativistic 



In the non-relativistic limit, (3(E) — y/2E/rrii. Fol- 
lowing the same procedure as above, we find that the 
non-relativistic form of each of —dx/dE and I is simply 
the first term of the corresponding relativistic relation. 



MIXED-SPECIES TARGETS 

For a target comprising multiple elements, the stop- 
ping power can be estimated from the combination of 
stopping powers from each element s 



dE y-^dE. , 
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- the Bragg addition rule. This relation is approximate 
because of chemical bond formation, which alters the ef- 
fective ionization. For ion energies much greater than the 
bond energies, the approximation should be accurate. 

Expanding as before and gathering terms, the range 
can be expressed very similarly to that for single-element 
targets. In the relativistic case, 
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where 



and 
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is the total electron density in the target. 
In the non-relativistic case, 
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which is again simply the first term of the relativistic 
relation. 



POLYNOMIAL FIT TO THE STOPPING 
DISTANCE SCALE 

Although well-characterized numerical approximations 
to the exponential integral exist [l7j |. they are not avail- 
able as standard functions in mainstream computer lan- 
guages, and require significant effort to implement from 
scratch. However, the logarithmic terms in the stopping 
power and range vary slowly compared with the powers 
of E; over wide ranges of energy, the stopping power can 
be approximated accurately by low-order polynomials. 
The use of polynomial approximations avoids the need 
to evaluate the exponential integral function. We define 
a stopping distance scale 



D = -Edx/dE, 



(15) 



which is particularly well-behaved above the low-energy 
singularity as it tends to zero with E, and increases 
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FIG. 1: Example calculation of stopping distance scale: pro- 
tons in water. The curve can be fitted well by a quadratic. 



monotonically (Fig. [T]). This quantity can be used as 
a crude, 0(1), (over)estimate of particle range, without 
requiring any series expansion or integration. Approxi- 
mating D by a polynomial 
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the range (Eq. [2]) is simply 
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ao must be zero, since E n / \n(aE) — » as E — > 0. 

To find the charged particle range in a specific sub- 
stance, it is straightforward to tabulate D(E) using Eq.[3] 
(and Eq. [10] for a multi-species target), fit a polynomial 
D P (E), and evaluate Eq. [T71 However, it is also possible 
to find a universal polynomial fit. Defining for conve- 
nience a different scaled energy 
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the stopping distance scale is (using Eq. [S]) 
D(F) 
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where the first term is the non-relativistic approxima- 
tion. The prefactor comprises universal constants and a 
simple problem-specific factor mj 2 j z 2 N Z . The relative 
magnitude of the relativistic term to the non-relativistic 
term depends only on /. Thus, by finding polynomial 
approximations 
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TABLE I: Polynomial fits to stopping distance scale functions. 



parameter 


5 < F < 100 


100 < F < 10000 


m 


2.17423 


1.64377 x 10 




2.29035 x 10 _1 


1.31696 x 10 _1 


n 3 


-3.36317 x 10" 4 


-4.55336 x 10~ 6 


n4 




2.07676 x 10" 10 




10 < F < 200 


500 < F < 10000 


T2 


-1.33946 x 10 


-9.92931 x 10 a 




2.04195 


3.49521 x 10 


re 


1.58588 x 10" 1 


1.00674 x 10 _1 


ra 


-7.67337 x 10~ 5 


-7.28447 x 10" 7 



(21) 



it is straightforward to find the polynomial coefficients 
for any problem: 
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and hence the range through Eq. [T7] 

The derivation presented above is valid for any choice 
of units. The Bethe relation breaks down when the log- 
arithm changes sign, i.e. when E approaches Irrii/4m e . 
Physically meaningful distances are obtained for greater 
energies. Using the Bloch estimate 12| for the effective 
ionization of the material, 



I ~ lOZq, 



(23) 



the relations are valid for E > 4600Zg. Polynomial fits 
were calculated over a range suitable for hadrons of en- 
ergy ~MeV to GeV (Table |J) . The relativistic term di- 
verges rapidly outside the fitting region. 

ENERGY DEPOSITION PROFILE 

Given E(x), the profile of energy deposition 
—dE{x)/dx (Bragg curve) can be calculated. E(x) is 
the inverse of the range, 1(E). 

For the non-relativistic range relation, E(x) can be ex- 
pressed in terms of the inverse of the exponential integral 
function, 



dE{x) q A z 2 NZ 
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If the stopping distance scale is represented locally in 
energy by a sufficiently simple polynomial, then it too 
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FIG. 2: Example calculation of Bragg curve via a polyno- 
mial (quadratic) fit to D(E): protons in water. The curve 
is presented backward from its usual form: as if the particles 
are accelerating from rest. Conceptually, the curve can be 
continued to arbitrarily high energies, i.e. long ranges. 

may be used to calculate —dE/dx. For example, taking 
a local quadratic fit 

D, p = axE + a 2 E 2 l p = a x E + y# 2 , (26) 

one obtains 

dE ~ * (27) 

dx y/aj + 2a 2 x 

(Fig.©. 

EXAMPLE CALCULATIONS 

Reference calculations are used to validate radiation 
protection simulations using different computer pro- 
grams. Here we compare the analytic solutions of the 
Bethe relation with results from widely-used programs 
SRIM [18], which uses numerical solutions of more de- 
tailed stopping powers developed from the Bethe rela- 
tion, and MCNP , which collects Monte-Carlo statis- 
tics for the simulated interaction of individual particles. 
Trial calculations were made for protons and a-particles 
stopping in Al and water. 

We use the Bloch estimate, Eq. [23] for the effective 
ionization of the target material. More accurate calcula- 
tions have been developed more recently, but the original 
Bloch estimate serves to demonstrate the correctness of 
our analysis. 

The results are consistent with the accuracy of the 
Bethe relation itself (Table HI)) , and are consistent with 
direct numerical integration of the Bethe equation with- 
out being affected by the low energy singularity. The 



TABLE II: Comparison between analytic calculation and 
computer simulations of ion ranges (in millimeters). 





anal vMf* 


MCNP 5 


SRIM 


p —> Al 








iu ivie v 


U.oy 


U.uz 


U.Do 


100 MeV 


36 


37 


37 


1 GeV 


180 


1510 


1530 


a — > Al 








10 MeV 


0.054 


0.062 


0.061 


100 MeV 


3.0 


3.2 


3.1 


1 GeV 


170 


180 


180 


p — > water 








10 MeV 


1.1 


1.2 


1.2 


100 MeV 


73 


77 


76 


1 GeV 


300 


320 


320 


a — > water 








10 MeV 


0.097 


0.110 


0.110 


100 MeV 


6.0 


6.3 


6.2 


1 GeV 


350 


380 


375 



greatest difference was for relativistic protons in Al. In 
this regime, radiative losses and nuclear reactions become 
significant [161 ] and the Bethe relation requires additional 
corrections. 

CONCLUSIONS 

Analytic solutions were found to power series expan- 
sions and polynomial fits to the Bethe relation. These so- 
lutions provide a convenient way to calculate ion ranges 
and energy deposition in regimes where the Bethe rela- 
tion is valid, i.e. kinetic energies of roughly 1-100 MeV/u, 
without depending on numerical integration. The use of 
a Taylor series restricts the accuracy at high energy; the 
relativistic expansion thus incorporates relativistic con- 
tributions to the range but is not valid to arbitrarily high 
energies. However, the analytic solutions can readily be 
used with more accurate formulations of the effective ion- 
ization. The accuracy was demonstrated by comparison 
with simulations from widely-used computer programs of 
ion ranges in Al and water. 
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